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Chiral effective field theory (EFT) enables a systematic description of low-energy hadronic inter¬ 
actions with controlled theoretical uncertainties. For strongly interacting systems, quantum Monte 
Carlo (QMC) methods provide some of the most accurate solutions, but they require as input local 
potentials. We have recently constructed local chiral nucleon-nucleon (NN) interactions up to next- 
to-next-to-leading order (N^LO). Chiral EFT naturally predicts consistent many-body forces. In 
this paper, we consider the leading chiral three-nucleon (3N) interactions in local form. These are 
included in auxiliary field diffusion Monte Carlo (AFDMC) simulations. We present results for the 
equation of state of neutron matter and for the energies and radii of neutron drops. In particular, 
we study the regulator dependence at the Hartree-Fock level and in AFDMC and find that present 
local regulators lead to less repulsion from 3N forces compared to the usual nonlocal regulators. 

PACS numbers: 21.60.Ka, 21.30.x, 21.65. Cd, 26.60.C 


I. INTRODUCTION 

Ghiral effective field theory (EFT) provides a system¬ 
atic expansion for nuclear forces based on the symmetries 
of QCD [iHS]. At a given order in the power counting, 
nuclear forces include contributions from pion exchanges 
and from shorter-range interactions. Ghiral EFT enables 
calculations with controlled theoretical uncertainties, a 
consistent description of electroweak interactions, and 
the matching to lattice QGD. In addition to nucleon- 
nucleon (NN) interactions, which are the dominant con¬ 
tribution to nuclear forces, chiral EFT naturally pre¬ 
dicts consistent many-body interactions, where the lead¬ 
ing three-nucleon (3N) forces enter at next-to-next-to- 
leading order (N^LO) [HIS]. It has been shown that 3N 
forces are key for the properties of neutron and nuclear 
matter [SHIll- A better understanding of 3N forces is a 
major frontier in nuclear physics. 

In addition to systematic chiral EFT interactions, re¬ 
liable many-body methods are needed. For strongly in¬ 
teracting systems, quantum Monte Carlo (QMC) meth¬ 
ods provide some of the most accurate solutions [m 
ng. These include Green’s function Monte Carlo 
(GFMC) [T7HT9] and auxiliary field diffusion Monte 
Carlo (AFDMC) [50] methods; for a recent review see 
Ref. [5T]. In continuum QMC calculations, the central 
object is the many-body propagator, which is of the form 

G(R, R'; St) = (R|e-^^-^|R') . (1) 

Here, R = (ri, si, r 2 , S 2 , • ■ •, tat, sn) is the configuration 
vector of all N particles, including the single-particle co¬ 
ordinates and spins r^, Si (and other quantum numbers). 
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St is a step in the imaginary-time evolution, and H is 
the Hamiltonian. 

In nuclear GFMC calculations all possible spin-isospin 
nucleon states need to be explicitly accounted for, which 
makes this method unsuitable for accurate neutron mat¬ 
ter studies due to an unfavorable scaling behavior. In 
contrast to GFMC, AFDMC rewrites the propagator by 
applying a Hubbard-Stratonovich transformation using 
auxiliary fields, which changes the scaling behavior fa¬ 
vorably at the cost of additional integrations over auxil¬ 
iary fields. We thus make use of the AFDMC method to 
study neutron matter. 

The trial wave function ipp in AFDMC is usually cho¬ 
sen to be of the form 


V^t(R) = 




i<j 


‘J>a(R), 


( 2 ) 


where inter-particle correlations are included through 
the Jastrow factor fj{rij) and is the noninteracting 
ground state given by a Slater determinant 


$a(R) = A 


4^ai (^ 2 , t^i) 


(3) 


where labels single-particle states (plane waves for 
neutron matter and Hartree-Fock orbitals for neutron 
drops [55]). 

For the evaluation of the propagator, it is necessary 
to be able to separate all momentum dependences up to 
quadratic terms. This can be done for local interactions, 
where the propagator for the momentum-dependent part 
is a Gaussian integral that can be evaluated analytically, 
while the interaction part can be easily obtained from 
the configuration vector (for more details see Ref. [21]). 
Chiral EFT interactions are naturally formulated in mo¬ 
mentum space and usually contain several sources of non¬ 
locality. 
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FIG. 1. Contributions to 3N forces at N^LO. These include a 
two-pion-exchange part given by the couplings ci, cs, and C 4 , 
a one-pion-exchange-contact interaction given by cd, and a 
3N contact interaction given by ce- 


Recently, local chiral NN potentials have been con¬ 
structed up to N^LO and have been used to calculate 
the energy of neutron matter and light nuclei using con¬ 
tinuum QMC methods P51 - E5] . Following the same strat¬ 
egy, a minimally nonlocal NN potential was developed in 
Ref. [26] with explicit A degrees of freedom. Monte Carlo 
methods have also been used to study neutron matter 
based on lattice techniques and other momentum-space 
QMC approaches [27H29|- 

For a complete calculation at N^LO, it is necessary to 
include 3N forces. In Sec. |Tl| we present local chiral 3N 
forces at N^LO, which are consistent with the local NN 
interactions of Refs. |23l|24]. The general expressions for 
the local 3N forces are given in the Appendix. We study 
in detail the regulator dependence of the leading two- 
pion-exchange 3N energy contributions at the Hartree- 
Fock level in Sec. IIIII and in AFDMC in Sec. IIVI This 
shows that present local regulators lead to less repulsion 
from 3N forces compared to using the usual nonlocal reg¬ 
ulators. We present results for the equation of state of 
neutron matter in Sec. |IV| a nd for the energies and radii 
of neutron drops in Sec. |V| Finally, we summarize and 
give an outlook. 


II. CHIRAL 3N FORCES IN COORDINATE 
SPACE 

In chiral EFT, the leading 3N forces at N^LO have 
three contributions: a two-pion-exchange part given by 
the couplings ci,C 3 , and C 4 , a one-pion-exchange-contact 
interaction given by cd, and a 3N contact interaction 
given by ce We show the N^LO 3N contributions 

diagrammatically in Fig. The Ci couplings are deter¬ 
mined by pion-nucleon or NN scattering, while Cd and 
Ce have to be fit to properties of A > 2 systems. 

Because we want to include 3N forces in AFDMC 
calculations in coordinate space, we have to use local 
coordinate-space expressions of the N^LO 3N forces, as 
was similarly done in Refs. |30l|3l]. To achieve this, we 
first Fourier transform the momentum-space expressions 
of the N^LO 3N forces. We begin with the 3N contact 
interaction Ve- In momentum space, this contribution 
vanishes in neutron matter due to the Pauli principle, 
when a regulator that is symmetric in the particle labels 
is used; see Ref. jSj. Because a local regulator does not 


fulfill this requirement, the Ve term will contribute. In 
this case, the regulator induces a finite range that mixes 
3N partial waves. After Fourier transformation, we find 
in neutron matter (with = 1) 

= (4) 


where we sum over all permutations 7r(ijfc) of the three 
particles i,j, and k, — r^, = 92.4 MeV is the 

pion decay constant, and we use A^ = 700 MeV. The 
expressions for general isospin and details on the Fourier 
transformation are provided in the Appendix. 

As for the Ve term, the one-pion-exchange-contact in¬ 
teraction Vd vanishes in momentum space for neutron 
matter due to the spin-isospin structure, if a symmetric 
regulator is used |6|. In coordinate space, the Vd term 
also contributes and after Fourier transformation we have 
two parts (see the Appendix): 


-rDjk _ 9a 


E 


Tr(ijk) 


47r 


S{r,j)X,k{rkj) 


- (Ti ■ crkS{rtj)6{rkj) , 


(5) 


where gA = 1.267 is the axial coupling, = 138.03 MeV 
is the averaged pion mass, and the function Xik{r) is 
given by 

X^k{r) = [S,kir)T{r) + a, ■ (Tk]Y{r ), ( 6 ) 

with the tensor operator S'ife(r) = 3<Ti-rcrfc-r—cTi-er*,, the 
function T(r) = l-|-3/(m^r)-|-3/(TO7rr)^, and the Yukawa 
function Y(r) = Here, r is the unit vector in 

the direction of r and r is the magnitude. There are two 
contributions from Vd, because one-pion exchange con¬ 
tains a long-range part as well as a delta-function part. 
The latter needs to be included to maintain the Gold- 
stone boson nature of the pion. 

We emphasize that there is an ambiguity in performing 
the Fourier transformation for the Vd term, depending on 
the choice of the initial spin-isospin structure. This leads 
either to terms S{rij)Xik{rkj) or 6{rij)Xkj{rkj) with dif¬ 
ferent spin indices in the X function. The two expres¬ 
sions are the same due to the <5 function, but lead to 
different results after regularization. Therefore, the dif¬ 
ferences from choosing different structures are a regula¬ 
tor effect and should be of higher order, as mentioned in 
Ref. [31] . These differences will vanish in the limit of in¬ 
finite momentum cutoff. Because we will not include the 
Vd term in our calculations in this paper, this effect will 
not influence the results here. The different Vd terms are 
studied in Ref. [5^ . 

We now turn to the two-pion-exchange contributions 
Vc in neutron matter (see the Appendix). For the part 
proportional to ci, we find 




ijk 


C.ci 


2/4(4x)2 


^ ^ Xij CTk ' Ykj 
■K{ijk) 
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xU{r.,)Y{n,)Uirk,)Y{rkj), (7) 


with the function U{r) = 1 + l/(m^r). This contribu¬ 
tion is similar to the long-range (LR) iS-wave part of the 
Illinois 3N forces; see Ref. [33]. 

The part proportional to C 3 is given by 


ijk _ 


TrlJK _ 

'^C,C3 - 


36/4 


E 

7r{ijk) 


4 2 

Tfl TTi 

X,j{r,j)Xkjirkj) - Xikirij)Sirkj) 


ml 


- -r^ Xik{rkj)5{v^j) Y a, ■ (Tk6{r,j) S{rkj) 
47r 


( 8 ) 


Thus, in coordinate space, for the C 3 part of the Vc 
term, there are four contributions due to a long-range 
and short-range part in each pion exchange. The 
first term ~ Xij(vij)Xkj{vkj) is a long-range two-pion- 
exchange contribution similar to the anticommutator 
part of the P-wave two-pion-exchange interaction of 
Ref. [33]. In addition, there is also a short-range (SR) 
term ~ 5{rij)5{rkj) which is similar to Ve but spin- 
dependent, and two intermediate-range (IR) terms ^ 
Xik{vij)5{rkj) + Xik{vkj)5{rij) similar to Vd- Note that 
although the spin-isospin structure is similar to the Ur- 
bana IX force and in general to the two-pion-exchange 
part of the Illinois forces, the spatial functions are quite 
different. Finally, the coordinate-space expression for the 
C 4 part of Vc is given in the Appendix. This does not 
contribute in neutron matter for general regulators due 
to the isospin structure. 

For a many-body system, the total 3N interactions are 
then given by Rsn = J2i<j<k Kn’ i,j,k = l,...,A. 
Moreover, in the AFDMC calculation is rewritten 
as a sum over cyclic permutations only. 

In order to regularize the local 3N forces consistently 
with the NN forces of Refs. [231 El]) we replace the <5 
functions by smeared-out delta functions of the form 




4 ^-{r/RsY^ 

7rF(3/4)Pi^ 


(9) 


where P 3 N is the three-body cutoff. For the long-range 
pion-exchange contributions, we multiply the Yukawa 
functions with the long-range regulator /long of Refs. [331 
[33] , given by 

r(r) ^ r(r)(^l - . ( 10 ) 


To be consistent with the NN cutoff Rq = 1.0 — 1.2 fm, 
we will also vary the 3N cutoff in this range, P 3 N = 1.0 — 
1.2 fm. We have checked that the IR and SR parts of Vc 
as well as the Ve and Ve> contributions in neutron matter 
vanish for P 3 N —> 0 (for infinite momentum cutoffs), and 
are therefore regulator effects. 



FIG. 2. Contributions to the neutron-matter energy per par¬ 
ticle E/N as a function of density n at the Hartree-Fock level. 
The black band shows the energy obtained using a nonlocal 
regulator, as in Ref. [34], with a 3N cutoff 400 — 500 MeV. 
The blue band corresponds to the LR part of the two-pion- 
exchange interaction Vc with the local regulator used here, 
the red band to the SR part of Vc, and the green band to the 
IR of Vc- For these bands, the cutoff in the local regulator 
is varied with R 3 N = 1.0 — 1.2 fm. The dashed-dotted line 
corresponds to the results for Vc using the local momentum- 
space regulator of Ref. m with a cutoff A 3 N = 500 MeV. 
This shows that these local 3N forces provide less repulsion 
at the Hartree-Fock level than with nonlocal regulators. The 
dashed lines show the results for Vc with the local regulator 
and small Rsn = 0.5 fm. 


In the following, we include all terms of Vc (the ci and 
C 3 parts for neutron matter), with the Ci couplings having 
the same values as in the local NN interactions [331 El] ■ 
Results including the shorter-range contributions Vd and 
Ve, which require fits of cd and ce, are studied in 
Ref. [33]. 


III. HARTREE-FOCK CALCULATION FOR 
NEUTRON MATTER 


We calculate the 3N contributions from the Vc part 
to neutron matter first at the Hartree-Fock (HF) level. 
This includes all interactions of Eqs. Q and Details 
on the HF calculation can be found in Refs. [34] . 

In Fig. [^ we show the contributions to the neutron- 
matter energy per particle E/N as a function of den¬ 
sity n. The blue band corresponds to the LR part of the 
two-pion-exchange interaction Vc with the local regula¬ 
tor used here, the red band to the SR part of Vc, and the 
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Rsn [fm] 

FIG. 3. Contributions to the energy per particle E/N at sat¬ 
uration density as a function of the cutoff i? 3 N. The lines 
show the LR, SR, and IR parts of the two-pion-exchange in¬ 
teraction Vc with the local regulator used here, calculated 
at the Hartree-Fock level. The bands are the contributions 
of the corresponding 3N parts to the AFDMC energies for a 
variation of the NN cutoff Ro = 1.0 — 1.2 fm; see also Fig. 


green band to the IR of Vc- For these bands, the cutoff in 
the local regulator is varied between i? 3 N = 1.0 — 1.2 fm. 
The dashed lines show the results for Vc with the local 
regulator and i? 3 N = 0.5 fm. In addition, the black band 
shows the energy obtained using a nonlocal regulator, as 
in Ref. [34], with a cutoff 400 — 500 MeV. 

The HF energy in neutron matter for the local Vc are 
in total « 3MeV at saturation density no = 0.16fm“^. 
This is only about half of the Vc contribution using the 
nonlocal regulator. The shorter-range (IR and SR) con¬ 
tributions, which are regulator effects, are small and with 
opposite sign. If we lower the coordinate-space cutoff, 
Rm = 0.5 fm (dashed lines), we find that the IR and 
SR parts almost vanish, as expected, and that the total 
HF energy is 5.5 MeV for the local Vc, which agrees well 
with the momentum-space result. We also note that the 
momentum-space result is very close to the infinite-cutoff 
result at the HF level. Thus, the smaller 3N energies for 
the local 3N forces are due to the local regulators used. 

To check this, we have performed a HF calculation of 
Vc using the local momentum-space regulator of Ref. |3T] 
with a cutoff of A 3 N = 500 MeV. This is given by the 
dashed-dotted line in Fig. At saturation density, we 
find an energy per particle of 3.8 MeV, which is com¬ 
parable to the result for the local 3N forces used here. 
This supports the above conclusion that 3N forces with 
local regulators provide less repulsion at the HF level 



FIG. 4. Variation of the AFDMG energy per particle at sat- 
nration density as a function of the 3N cutoff Rsn for an NN 
cutoff 1.0 fm (black lines in the upper part) and 1.2 fm (red 
lines in the lower part). The horizontal lines correspond to 
the NN-only energy. The squares are for the Ci and LR C 3 
part of Vc, the crosses include also the SR C 3 part of Vc, and 
the circles include all parts of Vc- 


compared to the usual nonlocal regulators. It may be 
possible that the Vjj and Ve parts, which contribute to 
neutron matter for local regulators, make up part of these 
differences. This is explored further in Ref. |32|. 

We show the Vc contributions to the energy per par¬ 
ticle E/N at saturation density as a function of the 3N 
cutoff i? 3 N in Fig. The lines show the LR, SR, and 
IR parts for the local regulator used here, calculated at 
the HF level. For all 3N cutoffs, the SR and IR parts 
are small and of opposite sign, while the major contri¬ 
bution of Vc comes from the LR parts. The SR and IR 
parts vanish for small coordinate-space (high momentum- 
space) cutoffs, as expected. The LR part increases up to 
the infinite-cutoff result at the HF level. In the cutoff 
range i? 3 N = 1.0 — 1.2 fm, the total contribution of Vc is 
« 3 MeV, and thus only about half of infinite-cutoff re¬ 
sult. This is what we also found in Fig.[^ We emphasize 
that the cutoff dependence from 400 — 500 MeV to infinte 
momentum-space cutoff is small for nonlocal regulators 
and these densities. 


IV. QMC CALCULATION FOR NEUTRON 
MATTER 

Next, we investigate 3N forces in neutron matter us¬ 
ing the AFDMC method, similarly to Refs. |2S1 [23] . In 
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FIG. 5. Dependence of the AFDMC energy per particle at 
saturation density as a function of the 3N cutoff i? 3 N on dif¬ 
ferent long-range regulators. Results are shown for an NN 
cutoff Rf) — 1.2 fm. The long-range regulator is given by 
[l — different parameters ni and n 2 . 


Fig.§ in addition to the HF results, we show the contri¬ 
butions of the LR, SR, and IR parts of Vc to the AFDMC 
energy, where the bands are from varying the NN cutoff 
i?o = 1-0 — 1.2 fm. For the SR part, the AFDMC en¬ 
ergies agree well with the HF energies, so that HF is a 
good approximation for this contribution. For the IR 
part, for large coordinate-space cutoffs, the agreement 
between HF and AFDMC results is good but worsens for 
smaller cutoffs. The uncertainty (from NN cutoff varia¬ 
tion) grows and the energy decreases significantly com¬ 
pared to the HF result. For the LR parts, the AFDMC 
energies are about 70 — 80% of the HF energies for higher 
cutoffs, which suggests that the LR N^LO 3N contribu¬ 
tions beyond HF are important. When lowering the 3N 
cutoff, the energy increases up to a plateau. Further 
lowering the cutoff, the system collapses and the energy 
rapidly decreases. In addition, the uncertainty grows. 

To study these effects more clearly, we show the vari¬ 
ation of the total AFDMC energy per particle at sat¬ 
uration density as a function of the 3N cutoff i? 3 N for 
an NN cutoff 1.0 fm (black lines in the upper part) and 
1.2 fm (red lines in the lower part) in Fig. The hor¬ 
izontal lines correspond to the NN-only energies. The 
squares include only the LR ci and C 3 part of He, the 
crosses include also the SR C 3 part of Ve, and the cir¬ 
cles include all parts of Vc- For the soft NN potential 
(i?o = 1.2 fm), we find the plateau of the AFDMC en¬ 
ergy to be at i? 3 N = 1.2 — 1.4 fm. If the cutoff is lowered, 
the energy decreases and for i? 3 N = 0.6 fm we find an 


FIG. 6. Energy per particle as a function of density for neu¬ 
tron matter at N^LO, including NN forces and the 3N Vc 
interaction in AFDMG. Results are shown for an NN cutoff 
Rq = 1.0 — 1.2 fm and 7?3 n in the same range. 

attractive 3N contribution^ For the harder NN poten¬ 
tial {Rq = 1.0fm), the plateau is found for smaller 3N 
cutoffs, i? 3 N = 1.0 — 1.2 fm. 

In general, the plateau is reached when ~ Rq, and 
the system collapses when i? 3 N is significantly smaller 
than Rq. This can be understood because harder NN 
potentials do not favor particles to be close and smaller 
3N cutoffs are needed to overcome this repulsion. If we 
want to decrease i? 3 N in our calculations, we also need 
to decrease Rq- Therefore, i? 3 N has to be chosen consis¬ 
tently with Rq which justifies our cutoff range. Because 
Rq < 1.0 fm is difficult for local NN potentials [24], we 
also do not decrease the 3N cutoff below that limit and 
use the range for Rq and i? 3 N within 1.0 — 1.2 fm. 

Although the collapse of the system for lower 3N cut¬ 
offs does not appear in the HF calculation, it is not an 
artifact of the AFDMC method. It is due to the function 
Xij{r), which includes terms ~ 1/r^. If three particles 
are in a small volume, this becomes very attractive, un¬ 
less it gets regularized with a large enough R 3 N. These 
cutoff values correspond to the position of the plateau. 
In a HF calculation, which only includes low-momentum 
states, this collapse will not appear. 

We have investigated the AFDMC energies when 


^ This behavior is qualitatively similar to the overbinding given 
by the Illinois 3N forces in pure neutron systems m- It would 
be interesting to see if using a similar cutoff would avoid the 
overbinding of neutron matter using Illinois forces. 
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FIG. 7. Comparison of the neutron-matter energy at N^LO 
based on the local chiral NN-I-3N potentials in AFDMC (this 
work) with the N^LO calculation of Ref. [10] based on the 
EGM N^LO potentials and using MBPT, with the particle- 
particle (pp) ladder results of Ref. |14| based on the EM N^LO 
potential, and with results based on the N^LOopt potential 
using CC theory and the SCGE method [Si- 


choosing different parameters in the long-range regula¬ 
tor function. In Fig. we show the dependence of the 
AFDMC energy per particle at saturation density as a 
function of the 3N cutoff i? 3 N on different long-range reg¬ 
ulators. Results are shown for an NN cutoff Rq = 1.2 fm. 
The long-range regulator is given by (l — 
with different parameters ni and n 2 - We find that the 
general picture is independent of the choice of the expo¬ 
nents in the regulator function. A consistent change of 
the short-range regulator has only a negligible effect on 
the energy. For different functions, the position of the 
plateau varies between 0.8 — 1.2 fm but the overall en¬ 
ergies at the plateau are comparable, generally ranging 
between 12.3 - 12.5 MeV. 

In Fig. we present the final result of our AFDMC 
simulations for the equation of state of neutron matter 
at N^LO. We show the energy per particle as a function of 
density including NN forces and the 3N Vc interaction. 
Results are shown for an NN cutoff Rq = 1.0 — 1.2 fm 
and i? 3 N in the same range. For the softer NN potential 
(i?o = 1-2 fm, lower lines) we find the energy per particle 
to be 12.3 — 12.5 MeV at saturation density for different 
3N cutoffs. The NN-only energy is 11.4 MeV and the 
3N Vc has an impact of « 1 MeV. For the harder NN 
potential (i?o = 1.0 fm, upper lines) we find an energy per 
particle of 15.5 — 15.6 MeV compared to 14.1 MeV for an 
NN-only calculation. Here, the impact of the 3N Vc is 


« 1.5 MeV. The variation of the total energy with the 3N 
cutoff is « 0.2 MeV in our cutoff range and considerably 
smaller than the variation with the NN cutoff, because 
the i? 3 N range lies in the plateau described above. 

We find the magnitude of the local 3N two-pion- 
exchange Vc forces to be at most about 1.5 MeV at sat¬ 
uration density, which is is smaller than a typical con¬ 
tribution of 4 MeV [Bj in momentum space with nonlocal 
regulators, including second- and third-order corrections. 
As discussed above, this difference can already be seen 
on the HF level and is most likely due to the present local 
regulators. This was also observed in the coupled-cluster 
calculations of Ref. m where a difference of 2 MeV was 
found for the neutron-matter energy per particle when 
choosing local versus nonlocal regulators with a similar 
cutoff of 500 MeV. Following these findings, local versus 
nonlocal regulators need to be further investigated. 

In Fig. we compare the neutron-matter energy 
at N^LO based on the local chiral NN-I-3N potentials 
in AFDMC (this work) with the N^LO calculation of 
Ref. (TD] based on the EGM N^LO potentials of Ref. [5B] 
and using many-body perturbation theory (MBPT), with 
the particle-particle (pp) ladder results of Ref. [14] based 
on the EM N^LO potential of Ref. [37], and with results 
based on the N^LOopt potential of Ref. [35] using the 
self-consistent Green’s function (SCGF) method [T3] and 
using coupled-cluster (CC) theory [TT]. At saturation 
density, the AFDMC energies are in general smaller than 
the other results, mainly due to the smaller contributions 
from local 3N forces. Furthermore, the density depen¬ 
dence of the AFDMC band is flatter than for the other 
calculations, which may be explained by differences in the 
NN phase shift predictions. We would expect the results 
to come closer when including chiral forces at next-to- 
next-to-next-to-leading order (N^LO). A comparison of 
AFDMC results with MBPT results using the same lo¬ 
cal potential, as in Refs. [331 (33], be presented in a 
forthcoming paper. 


V. NEUTRON DROPS 

Neutron drops in external potentials provide useful 
constraints for energy-density functionals and their ap¬ 
plications to neutron-rich nuclei [3^ [31] . They constitute 
a simplified model of neutron-rich nuclei, where the ex¬ 
ternal well simulates the effects of the core on the valence 
neutrons. Their study is therefore a natural addition to 
homogenous neutron matter. 

We have performed AFDMC calculations for the ener¬ 
gies and radii of neutron drops with TV = 8, 20, 40, and 
70 neutrons in a harmonic oscillator potential with an os¬ 
cillator parameter Hu = 10 MeV. For these calculations 
we used the same cutoff in the NN and 3N interactions, 
Rq = i? 3 N = 1.0—1.2 fm. The results for the energies and 
radii are tabulated in Table [I] and shown in Fig. We 
give the results at different orders in the chiral expansion, 
and at N^LO, for NN forces only, plus only the LR ci and 
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FIG. 8. Energies and radii of neutron drops with N = 8 and 20 neutrons in a harmonic oscillator potential with an oscillator 
parameter huj — lOMeV using AFDMC. The same cutoff is used in the NN and 3N interactions. We give the results at different 
orders in the chiral expansion, and at N^LO, for NN forces only, plus only the LR ci and ca parts of Vc, and including also the 
SR and IR parts of Vc- The bands are given by the cutoff variation Ro = Rsn = 1.0 — 1.2 fm. At LO with the softer cutoff and 
N = 20, the system collapses. We compare our results with the calculations of Ref. [40], using coupled-cluster theory at the 
ACCSD level, where the band is given by two different similarity renormalization group (SRG) evolution scales (for one initial 
NN-I-3N Hamiltonian). 


C 3 parts of Vc, and including also the SR and IR parts of 
Vc- The bands in Fig. are given by the cutoff variation 
Ro = Ron = 1.0 —1.2 fm. We generally find a good order- 
by-order convergence of the energies and radii. The band 
increases with neutron number; at the level of N^LO-I-Vc' 
it is 1% for the energy oi N = 8 neutron drops, 2% for 
N = 20, 5% for N = 40, and 7% for N = 70. Fur¬ 
thermore, in systems with N > 20 and at low orders, 
our calculations do not converge and a collapse occurs 
(for the interactions not listed in the table). This is due 
to the higher densities inside the larger neutron drops 
which are also tabulated for the N^LO-I-Vc' Hamiltonian 
at r = 0.125 fm in Table I. These show that the different 
particle numbers probe a broad range of central densities 
from low densities for N = 8 to twice nuclear saturation 
density for N = 70, connecting the neutron drop results 
with our neutron matter calculations. The higher the 
(central) density, the larger the effect of the 3N forces, 
leading to a collapse for high densities. 

Furthermore, for all N at N^LO the relative contribu¬ 
tion of Vc is always larger for the Rq = 1.2 fm potential 
than for the 1.0 fm one. For N = 70 we find Vc to con¬ 
tribute 3.7% for 1.2 fm versus 2.7% for 1.0 fm and for 
N = 8 0.70% versus 0.56%. This result is the opposite 
of what Fig. 6 shows for homogeneous matter: there the 
softer NN potential leads to a smaller 3N contribution. 
This is due to the higher central densities of the neu¬ 


tron drops for the softer potentials, leading to larger 3N 
contributions. 

In Fig. we also compare our results at N^LO with 
the calculations of Ref. HO] using coupled-cluster theory 
at the ACCSD level. The latter results are based on an 
SRG-evolved chiral Hamiltonian starting from the N^LO 
NN potential of Ref. [Hj with a cutoff of 500 MeV and 
local N^LO 3N forces, regulated in momentum space with 
the same cutoff value, including also the Vd and Ve parts. 
The band is given by two different SRG evolution scales 
(for this initial NN-I-3N Hamiltonian). We find very good 
agreement between the two approaches after inclusion of 
N^LO 3N forces, whose contribution is small. 


VI. SUMMARY AND OUTLOOK 

We have presented local chiral 3N forces at N^LO 
that are consistent with the local NN interactions of 
Refs. [231 mi- We have investigated the 3N two-pion- 
exchange contributions to neutron matter both at the 
HF level and in AFDMC calculations, including a de¬ 
tailed study of the regulator dependence. Our results 
show that present local regulators lead to less repulsion 
from 3N forces compared to using the usual nonlocal reg¬ 
ulators. This is already present at the HF level. 

In neutron matter, the dependence on the 3N cutoff 









TABLE I. Energies (in MeV) and radii (in fm) of neutron 
drops with N = 8, 20, 40, and 70 neutrons in a harmonic os¬ 
cillator potential with an oscillator parameter huj = 10 MeV. 
The same cutoff is used in the NN and 3N interactions. We 
give the results at different orders in the chiral expansion, and 
at N^LO, for NN forces only, plus only the LR ci and C3 parts 
of Vc, and including also the SR and IR parts of Vc- For the 
latter, we also give the central densities (at r = 0.125 fm in 
fm“®). In systems with V ^ 20 a collapse occurs at low orders 
(for the interactions not listed in the table). 


N 

Hamiltonian 

E 

rms radius 

ric 

8 

LO(l.O) 

133.51(3) 

2.571(2) 


8 

NLO(l.O) 

133.31(3) 

2.633(2) 


8 

N2LO(1.0) NN-only 

134.83(3) 

2.656(2) 


8 

N2LO(1.0)-tVc.LR 

135.53(3) 

2.673(2) 


8 

N2LO(1.0)-tV'c 

135.59(5) 

2.677(2) 

0.07(1) 

8 

LO(1.2) 

123.87(8) 

2.365(2) 


8 

NLO(1.2) 

133.11(3) 

2.607(2) 


8 

N2L0(1.2) NN-only 

133.86(3) 

2.616(2) 


8 

N2LO(1.2)-tV'c.LR 

134.78(3) 

2.638(2) 


8 

N2LO(1.2)-tV'c 

134.80(5) 

2.637(2) 

0.08(1) 

20 

LO(l.O) 

432.92(12) 

2.966(2) 


20 

NLO(l.O) 

428.67(8) 

3.062(2) 


20 

N2LO(1.0) NN-only 

434.93(6) 

3.089(2) 


20 

N2LO(1.0)-tV'c.LR 

440.73(7) 

3.137(2) 


20 

N2LO(1.0)-tV'c 

440.67(10) 

3.138(4) 

0.15(1) 

20 

NLO(1.2) 

422.13(5) 

2.987(2) 


20 

N2L0(1.2) NN-only 

424.64(5) 

2.989(2) 


20 

N2LO(1.2)-tV'c,LR 

432.05(5) 

3.041(2) 


20 

N2LO(1.2)-tUc 

431.55(10) 

3.036(2) 

0.18(1) 

40 

NLO(l.O) 

1056.24(26) 

3.459(2) 


40 

N2LO(1.0) NN-only 

1071.17(24) 

3.481(2) 


40 

N2LO(1.0)-tVc.LR 

1094.90(18) 

3.557(2) 


40 

N2LO(1.0)-tV'c 

1093.18(32) 

3.556(2) 

0.16(1) 

40 

NLO(1.2) 

1018.25(25) 

3.318(2) 


40 

N2L0(1.2) NN-only 

1019.44(20) 

3.293(2) 


40 

N'‘LO(1.2)-tVc,LR 

1048.70(20) 

3.385(2) 


40 

N2LO(1.2)-tUc 

1045.81(50) 

3.377(2) 

0.20(1) 

70 

N2LO(1.0) NN-only 

2235.89(60) 

3.877(2) 


70 

N2LO(1.0)-tV'c.LR 

2302.05(150) 

3.987(2) 


70 

N2LO(1.0)-tV'c 

2296.23(111) 

3.991(2) 

0.25(2) 

70 

N2L0(1.2) NN-only 

2071.90(70) 

3.593(2) 


70 

N2LO(1.2)-bVc.LR 

2164.53(180) 

3.730(4) 


70 

N2LO(1.2)-tUc 

2148.75(190) 

3.711(4) 

0.31(2) 


over the range i? 3 N = 1.0 — 1.2 fm is small compared 
to the NN cutoff variation, but for lower 3N cutoffs the 
system starts to collapse. We have also studied the influ¬ 
ence of different local 3N regulators and found that the 
general picture remains the same. Our findings lead to 
the conclusion that local versus nonlocal regulators have 
to be extensively studied. It will be crucial to develop a 


method of assessing the quality of local regulators to find 
improved versions for these regulators. 

We have studied the neutron-matter equation of state 
for local chiral NN and 3N interactions and found smaller 
energies compared to other calculations, mainly due to 
less repulsion from 3N forces. We also simulated neutron 
drops in an external harmonic oscillator potential for neu¬ 
tron number N = 8, 20, 40 and 70 and investigated their 
energies and radii for different chiral orders. Our results 
show very good agreement with previous coupled-cluster 
calculations using chiral potentials [JD] (also with local 
3N forces). 

Work on the determination of the two 3N couplings 
CD and ce for the local 3N forces is reported in Ref. [32] . 
The inclusion of the full N^LO 3N forces will enable novel 
many-body calculations of nuclei and nuclear matter with 
QMC methods based on chiral EFT interactions. 
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Appendix A: Coordinate-space expressions 


In momentum space, the N^LO 3N interactions are 
given by mil] 

2 


Uc = 


1 / 5T 


2 V2/. 




E 


(Al) 


QA CD 


8/2/2A^^E 


E • Qfc T* • Tfe , (A2) 


Ve = 


Ce 


2/^A^ 


^ ^ '^i ' '^k ■. 
■n-(ijk) 


(A3) 


where = Pi — Pi is the momentum transfer of parti¬ 
cle i (all other quantities are defined in Sec. and 
includes the different contributions from the Ci’s 


_ xa/3 
^ijk - 0 


4cito^ 


E C4 
J Jtt 


+ 


2c3 

J TT 


J2 ' f2 


T-j CTj 


(q* X qfc) ■ 


(A4) 
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In neutron matter, the 3N contributions simplify because 
the isospin structure can be evaluated explicitly, with all 
Ti ■ Tj = I and the C 4 part vanishes [5] . 

We Fourier transform the 3N interactions with respect 
to the momentum transfers of particles i and k, which 
yields the coordinate-space expression as a function 
of rij and Because the 3N interactions include a 
sum over all permutations, taking a different choice for 
the momentum transfers would lead to the same result. 
However, this will not be the case when a regulator in 
momentum space is included before Fourier transform¬ 
ing. For the Ve contribution this gives 


ytjk ^ d^qk 

J (27r)3(27r)3 

Ce 


Ve , 


= ^ d{vj)6{rkj ). (A5) 

For the Fourier transformation of the Vd contribution 
one has 


(27r)3 (27r)3 
CogA 


mK 




77{ijk) 


d^qi 

(27r)3 




X / (A 6 ) 


(27r) 


% 


The second integral gives an expression similar to one- 
pion exchange: 


(Pqu o-fe • Qfe I 


qfc 


(27r) 

o 2 


ql + ml 




7Tl^ 1 

= - TTrXikirkj) + - (T^ ■ (TkS{rkj), (A7) 

iZTT o 

with Aifc(r) defined in Eq. (§. As a result, in addition 
to the one-pion-exchange-contact part, the Fourier trans¬ 
formation also leads to a contact-contact part in 


= -iaTd^ 


m. 


47rr, 


= ^^a^f%U{n,)Y{n,). (AlO) 


This leads to 


^r^ok _ Cim^gA V- __~ ^ 

^C,ci - p /•4('4_\2 CTfc • Tfcj 

X U{n^)Yin^)Uirkj)Y{rkj). (All) 

Next, the Fourier transformation of the C 3 part of Vc 
gives 


yijk _ CsgA ^ ^ 

^C,c3 ~ A fi ^ 

^(ijk) 

(fqk (Tk ■ Qfc 


(Pqi (Ti ■ q, 
(27r)3 qf + ml 




(A12) 


J (27r)3q2 + TO2 

Similar to the Fourier transformation for one-pion ex¬ 


change in Eq. (A7) one gets 


d^qi (Ti ■ Qj 

(27r)3 qf + ml 




= -?^f 

47r 


’'rx«-3^“')Uni)r(r«) 




Combining this leads to 


yijk _ 

C.C3 - 36^4 


E 

7r{ijk) 


' '^k 


(A13) 


(A14) 


4 2 

TTL Tfl 

Xij (v^j)Xkj (vkj) - ^Xik {v^j)5 {vkj) 


_(47r) 


47r 


- -r^ Xik{Vkj)5{Vij) + (Ti ■ (Tk5{Yij)5{Vkj) 
47r 

Finally, for the C 4 part of Vc we obtain 


^ijk _ ^D9A 

^ - 24/4A, 


E 

TT{ljk) 


T* • Tk 


47r 


S{nj)Xik{rkj) 


- ■ crkd{rij)6{rkj) 


IJJ 

For the ci part of Vc we have 

.^ijk _ cimlg\ r d^q^ (Ti ■ 

2/4 ^2Tf,l^ml 

f d^qk (Tk-qk 
J {2Tr ql+ml 

The integrals are readily evaluated using 


(AS) 


(A9) 


d^qi (Ti ■ qj 

(27r)3 ql + ml 






Tr(ijk) 

^4 


2i(47r) 


^[Xijirij),Xkj{rkj)] 


- ■ Ak X (T )(1 - T{rC)yirij)5{rkj) 

47r 

9 

TTi^ 

- -^o■^ ■ Ak X (Tj) (1 - T{rkj))Y{rkj)5{vC 


3^ 

47r 

47r 


(Ti ■ r,n r, 


{(Tk X (Tj)T{r^J)Y{rij)5{Ykj) 


’’k ■ Ykj Ykj ■ {(Tj X (Ti)T{rkj)Y{rkj)S{rC 


+ cTt ■ {(^k X (Tj)S{rij)6{rkj) 


(A15) 
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